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Abstract. Ultra High Energy (UHE) particles coming from discrete extra- 
galactic sources are potential candidates for EAS events above a few tens of 
EeV. In particular, galaxies with huge infrared luminosity triggered by collision 
and merging processes are possible sites of UHECR acceleration. Here we check 
whether this could be the case. Using the PSCz catalogue of IR galaxies we 
calculate a large scale anisotropy of UHE protons originating in the population of 
the luminous infrared galaxies (LIRGs). Small angle particle scattering in weak 
irregular extragalactic magnetic fields as well as deflection by regular Galactic 
field are taken into account. We give analytical formulae for deflection angles with 
included energy losses on cosmic microwave background (CMB). The hypotheses 
of the anisotropic and isotropic distributions of the experimental data above 40 
EeV from AGASA are checked, using various statistical tests. The tests applied 
for the large scale data distribution are not conclusive in distinguishing between 
isotropy and our origin scenario for the available small data sample. However, 
we show that on the basis of the small scale clustering analysis there is a much 
better correlation of the UHECRs data below GZK cut-off with the predictions of 
the LIRG origin than with those of isotropy. We derive analytical formulae for a 
probability of a given number of doublets, triplets and quadruplets for any density 
distribution of independent events on the sky. The famous AGASA UHE triple 
event is found to be very well correlated on the sky with the brightest extragalactic 
infrared source within 70 Mpc - merger galaxies Arp 299 (NGC 3690 + IC 694). 



1. Introduction 



The existence of the UHE cosmic rays, after their discovery almost half a century 
ago, remains still puzzling. In advent of the new giant experiment named in honour 
of Pierre Auger, there are many proposals to explain their origin. UHE cosmic rays 
(UHECR) seem to be extragalactic because their arrival directions are not correlated 
with the Galactic disk. Thus, sites of their origin should be different from normal 
galaxies (like our Galaxy). In this work we check the hypothesis that the powerful 
luminous infrared galaxies (LIRGs) might be the UHECR sources. Large fraction of 
bright IR galaxies are found to be interacting systems, suggesting that collision and 
merging processes are mostly responsible in triggering the huge IR light emission . 
Fraction of interacting systems increases with IR luminosity and in the population 
of the most IR luminous objects in the Universe almost all appear as gravitationally 
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interacting. Favourable environments for accelerating particles to UHE regime via 
the first order Fermi process are provided by amplified magnetic fields on the scale of 
tens kpc resulting from gravitational compression, as well as high relative velocities 
of galaxies and/or superwinds from multiple supernovae explosions ^. There 
have been several observational claims that colliding galaxies could be possible sites 
of the UHECR origin. Al-Dargazelli et al [I have proposed that some clustering 
of UHECR shower directions above 10 EeV| are associated with nearby colliding 
galaxies. Uchihori et al have analysed combined world data from Northern 
hemisphere experiments and concluded that two triplets and a doublet lie in a vicinity 
of interacting systems. Takeda et al ||2^ have noticed that the interacting galaxy VV 
141 at z=0.02 is a possible candidate for the triplet of events above 50 EeV from the 
AGASA experiment. However, as we shall show, there is another favourable candidate 
for the origin of this UHECR cluster - Arp 299 (Mrkl71, VV 118a/b), a member of 
LIRG class of extragalactic objects. This system (RA= 171.4°, 6 — 58.8°), consisting 
of two interacting starburst galaxies, is the closest extragalactic object (distance 
42 Mpc for Hq = 75km s~^ Mpc~^) with IR luminosity greater than 5 x lO^'^^L© 
(~ 2 X W^^ergs s~^) and it is the brightest IR source within 70 Mpc. Such high 
IR luminosity is related to young and violent star forming regions. There is also 
observational evidence of superwind outflows, large scale strong radio emissions and 
the estimated supernova rate is about 0.6 per year |||, |l^, 14, 2^. 

The aim of the present paper is to check whether the LIRGs could be the sites 
of origin of the UHECRs observed at the Earth. Using the PSCz catalogue |20 we 



construct the all-sky maps of UHE proton intensities originating in LIRGs, taking into 
account effects of particle propagation through the extragalactic medium and, as an 
example, possible influence of the regular galactic magnetic field (GMF). We check 
both hypotheses: origin in LIRGs and, on the other hand, the isotropic distribution of 
the experimental data above 40 EeV from AGASA |n2[, using various statistical tests. 



2. Anisotropy calculations 

2.1. The PSCz catalogue 

The PSCz catalogue consists of almost 15000 IR galaxies with known redshifts, 
covering 84% of the sky. It should be noted that observational limitations cause that 
some of the extragalactic objects may not have measured redshifts or may be even 
unobserved in the dust obscured regions within the disk of the Galaxy (the so called 
zone of avoidance) . Figure 1 (top) presents the distribution of all PSCz galaxies with 
known redshifts. We can see patches on the sky regions excluded from the catalogue. 
On inspection of the superimposed directions of the AGASA showers above 40 EeV 
we can find only a few cases where experimental events lie in the vicinity of excluded 
regions, which should not affect much our analysis. As the potential UHECR sources 
we have selected objects with luminosities in the far infrared (FIR) range exceeding 
LIRG limit Lpjn = lO^^L©, that is over one order of magnitude greater than the 
estimated FIR luminosity of our Galaxy |^], and with distances up to 1 Gpc, giving 
the total number of 2811 sources (figure 1, bottom). A large fraction of such objects 
show to be in an apparent stage of collision and merging. 



t lEeV = 1 X IQi** eV 
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Figure 1. Top: Sky distribution of 14677 PSCz galaxies with known redshifts 
and superimposed AGASA showers above 40 EeV (stars). Aitoff projection in 
equatorial coordinates. Bottom: Selected 2811 objects with Lpjn > W^^Lq and 
with distances up to 1 Gpc 
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2.2. CR propagation 

To predict the CR anisotropy expected in the model of LIRG origin we have assumed 
the foUowing: 

• Protons are injected at the sources with a spectrum dN/dE ^ truncated at 
10^1 eV. 

• CR luminosity at the source is proportional to its L pjji luminosity. 

• CR propagate through the intergalactic medium, where they are scattered and 
suffer energy losses on the cosmic microwave background (CMB)g. 

To calculate CR arrival directions we have derived an analytical formula for the 
distribution of the deflection angles, for the multiple small angle scattering in 
weak irregular magnetic fields with continuous energy losses taken into account (see 
Appendix A). Since strength and structure of the extragalactic magnetic field are 
largely unknown, we adopt here the upper limit for this field magnitude and coherence 
scale: B (^c)^^^ = 1 nG (IMpc)^/^, as measured by Faraday rotation of radio signals 
from distant quasars ||l5|. Because the analysis is strongly energy dependent, we 
consider two energy regions 40 to 80 EeV and above 80 EeV, where this limit is just 
about the Greisen-Zatsepin-Kuzmin (GZK) flux cut-off predicted from CR interactions 
with the 2.7 K CMB radiation ||. 

2.3. Maps of the expected anisotropy 

Equatorial maps of the expected intensity of the UHE protons originating in LIRGs for 
two energy regions, 40-80 EeV and above 80 EeV, with sky coverage and declination 
dependent exposure for the AGASA experiment |^ (see also Appendix B), are 
presented in figures 2 (top) and 3 with superimposed AGASA events. In figure 2 (top) 
we can see that expected proton intensities for 40-80 EeV show good correlation with 
the distribution of the experimental events from AGASA. Especially, clearly visible is 
the region of the sky where high proton intensities from Arp 299 correlate with the 
AGASA triplet of events (RA« 170°, (5 « 60°). At higher energies, because of the 
rapid proton energy losses on CMB, only sources located at distances not larger than 
100 Mpc are visible on the map. No correlation of the data events above 80 EeV with 
the expected CR intensities can be seen (figure 3). There is an apparent group of 
UHE showers in the region with RA« 300°, where no possible sources exist. However 
some of the showers from this group lie close to the region heavily obscured by the 
Galactic disk. 

2.4. Influence of Galactic magnetic field 

The global regular GMF structure is not well known. Analyses of rotation measures 
of pulsars and extragalactic radio sources suggest that the GMF has a bisymmetric 
spiral (BSS) form with field directions reversed from arm to arm [^2| , |ll[. This is 
also supported by the observed value of the pitch angle of the local field and number 
of field reversals within the Galactic disk |ll|, |l^. To examine the influence of 
the regular GMF on the extragalactic UHE fluxes we have chosen the bisymmetric 
spiral model with field reversals and odd parity (BSS- A) adopted from Stanev p^ . 
Because the predictions are strongly model dependent, and the procedure applied 
gives a rather simplified picture, the analysis presented here is only an example of the 
possible extragalactic CR flux distortion. It is presented for the lower energy range 
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Figure 2. Expected maps of proton intensities without (top) and witii (bottom) 
influence of the regular galactic magnetic field (model BSS-A) for protons with 
40-80 EeV originating in LIRGs to be seen by AGASA, with superimposed 
47 AGASA shower directions in this energy range (stars scaled with energy). 
Contours of constant flux per unit solid angle are spaced linearly. 
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Figure 3. Expected map of proton intensities (energy above 80 EeV) from LIRGs 
to be seen by AGASA, with superimposed 11 AGASA showers in this energy range 
(no GMF). 

40-80 EeV, where the effect is stronger and easily visible. 

Simulations of a large number of monoenergetic antiprotons ejected from the 
Earth and followed to the halo border through GMF give flux modification factors 
depending on directions in the 'extragalactic' sky. Then, for an assumed flux 
distribution of protons at the Galactic borders we are able to calculate their fluxes 
at the Earth. Distortion effects i.e. reduction or magnification and shifting of the 
particle fluxes are visible in Figure 2 (bottom) . The modification of the extragalactic 
UHECR intensities due to the influence of the GMF BSS-A model seems to worsen 
the correlation with the experimental data seen in the upper figure. 



3. Statistical tests and results 
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3.1. Smirnov-Cramer-von Mises test 

To check the hypotheses of isotropic and anisotropic distribution of the experimental 
data above 40 EeV from AGASA, we have used Smirnov-Cramer-von Mises (SCvM) 
free of binning test modified for a 2-dimensional distribution analysis. SCvM test 
is based on comparing the cumulative distribution function F{X) under hypothesis 
Ho with the equivalent function Sn{X) of the data. The considered statistics is 
defined as follows 

W^^ (SN{X)-FiX)) fiX)dX (1) 

J —oo 

where f{X) is the probability density function corresponding to the hypothesis Hq. 
Sn{X) is based on the experimental data (two coordinates of the UHECR events on 
the equatorial map) and always increases in steps of equal height, N~^, where A'' is 
the total number of data. It is worth to note that this test is reliable even for small 
statistics {N > 3). Critical NW'^ values for the confidence level a —0.1 and 0.05 are 
0.347 and 0.461 respectively. 



E 


Tested hypotheses 


NW^ values 


40- 


Isotropy 


0.063; 0.116; 0.202; 0.263 


80 


LIRGs anisotropy 


0.069; 0.121; 0.106; 0.215 


EeV 


LIRGs anisotropy + GMF 


0.115; 0.195; 0.181; 0.247 


>80 
EeV 


Isotropy 

LIRGs anisotropy 


0.042; 0.573; 0.658; 0.071 
0.631; 0.351; 2.177; 0.338 



Table 1. 



Each of the cumulative distribution functions, F{X) and Sn{X), is the integral 
of the probability density function over the rectangle area defined by the coordinate 
point and one of the corners of the map. In this way, from the 2-dimensional maps 
we have constructed a 1-dimensional F{X) and Sn{X), for four cases depending 
on the chosen corner of the map. Although the four obtained NW"^ values are not 
completely independent, they provide a valuable insight into the analysis. Results of 
tested hypotheses: isotropy and LIRGs anisotropy (with and without GMF), energies 
40-80 EeV and above 80 EeV are shown in table 1. All the hypotheses in the energy 
range 40-80 EeV pass, assuming confidence level a =0.1 (A^W^^=0.347). However, 
with the presence of the GMF the agreement is worse. The hypotheses for energies 
above 80 EeV fail, especially in LIRG anisotropy scenario. 

Irrespective of the above, table 1 shows that it is not easy to draw conclusions 
from this test as the NW'^ values calculated for different map corners scatter quite 
significantly. Thus, it is advisable to apply another, hopefully more powerful, 
statistical test 

3.2. Eigenvector test 

The idea is based on assigning unit directional vectors to the data points on the 
celestial sphere. By finding the normalized eigenvalues ti, T2, T3 of the orientation 
matrix T constructed for the N unit data vectors it is possible to discriminate between 
the isotropic and some anisotropic distributions Assuming < ti < T2 < rs < 1 
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Figure 4. Left: Eigenvalue analysis distributions (each point corresponds to 
47 events) drawn from isotropy (dots) and LIRGs anisotropy (crosses) scenarios. 
Points based on AGASA data (star and thick cross represents respectively our 
and Medina Tanco calculation from 'old' AGASA data - 47 events with energy 
above 40 EeV. The open square is for 'new' data - 47 events, 40-80 EeV); Right: 
Simulations, as previously, but for large samples of 500 events each. 



the empirical shape criterion 7 = [logio{T:i/T2)] / [logio{T2/Ti)], and the strength 
parameter ^ = /ogio(T3/ri) are used to discriminate girdle type from clustered 
distribution. Distributions of the girdle and cluster type plot with 7 below and above 
unity, respectively. Intermediate i.e. partly girdle, partly cluster distributions plot 
around the line 7=1. Isotropic distributions plot with strength C near zero. 

In figure 4 (left) there are shown distributions of 500 samples, consisting of 
47 events each, simulated under isotropic (dots) and LIRGs anisotropic (crosses) 
scenarios. The open square denotes AGASA data, 47 events with energy 40-80 EeV. 
The two hypotheses, LIRG anisotropy and isotropy, correspond to the two regions 
in the C ~ 7 plane which overlap significantly for the 47 event samples. Mostly due 
to the small statistics considered, the eigenvector analysis is not sensitive enough to 
distinguish between isotropic and LIRG anisotropic origin of the experimental data. 
Hopefully with a sample of more than 500 events, to be collected in a few months of a 
full operational mode of the Auger observatory, the two distributions should separate 
(figure 4, right). 

The eigenvector method has been already applied to the UHECR anisotropy 
analysis by Medina Tanco The point corresponding to the 'old' AGASA data 

(47 events above 40 EeV), calculated by this author, is shown on the C ~ 7 plane 
(thick cross in figure 4, left) lying well away from the simulated isotropic distribution. 
However, the point calculated by us for the same AGASA data (star) has a lower 
strength value and lies on the isotropy distribution. Thus, contrary to the result of 
Medina Tanco, the AGASA data do not show such a strong clustering (on the basis 
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Figure 5. Probability distributions of doublets for isotropic (circles) and LIRGs 
anisotropic (crosses) scenarios for the cases of zero (solid line) and one triplet 
(dashed line). 

of the eigenvector analysis) to differ from the isotropy distribution. We have checked 
the correctness of our result by calculating the eigenvalues of the orientation matrix 
also analytically (which is possible for a 3x3 matrix). 

We have not done the eigenvector analysis for the data above 80 EeV because of 
the sample of 11 events was too small. 

3.3. Multiplet analysis 

Apart from analysing the large scale LIRG anisotropy (as in previous chapters) we 
will check this hypothesis by analysing probabilities of multiplets, i.e. groups of events 
with small angular separation. 

In the available|| AGASA data from 40-80 EeV there are two doublets (separation 
angles smaller than 2.5°) and one triplet (with criteria taken from the AGASA 
pubhcation H^, see also Appendix B). Figure 5 shows probability distribution of 
the number of doublets based on the 10* samples of 47 events each, simulated from 
the LIRGs anisotropic map (crosses) for the cases of zero and one triplet. There 
are also similar probabilities for the isotropy scenario (circles), calculated analytically 
(Appendix B). From figure 5 we find that the probability of obtaining at least the 
observed number of multiplets from LIRGs anisotropy scenario is 0.129 (7.1 x 10~^) 
for more than two (three) doublets, while from isotropy (with the same assumed 
response in declination of the AGASA experiment) it is only 0.67 x 10~^ (0.2 x 10^^). 
All AGASA events above 40 EeV have one triplet and six doublets, a collection even 

§ AGASA claims recently pa] three doublets in the energy range (4-10)10i^eV but it is not clear 
whether the energy of one shower is below 810^*' eV 
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Figure 6. Smoothed cos(9) distributions for separation angles between AGASA 
data (47 events 40-80 EeV) and 5000 events drawn from isotropic (dotted line) and 
LIRG anisotropic (solid line) distributions. Dashed line represents distribution of 
cosines within simulated LIRG events (47x5000) averaged over 10^ samplings. 

less probable to be obtained from isotropy. 

The above analysis may indicate only that there is some anisotropy (which 
increases mainly the probability of triplet) but not necessarily that correlated with 
LIRGs. A relationship between the set of the data and our hypothesis can be checked 
by looking for small angle correlation between the simulated and the real events. In 
figure 6 there are shown smoothed distributions of cos9 in the range 0.9—1, where ^'s 
are separation angles between the AGASA events, from one side (47 events, 40-80 EeV) 
and, from another side, 5000 events drawn from the isotropic (dotted line) or from 
the LIRG (solid line) distributions (giving 47x5000 cosine values for each case). The 
obvious flat distribution of cos9 between the AGASA data and isotropic events (dotted 
line) with a small rise from left to right shows the effect of the AGASA declination 
efficiency. The dashed line represents cos9 distribution for taken between the events 
simulated from the LIRG map (47x5000 events), averaged over 10-'' such runs. Here, 
a significant rise in the range of a few degrees {9 < 8 deg, cos9 > 0.99) indicates a 
clustering of the simulated LIRG events. A similar rise is seen in the distribution of 
cos9 between AGASA data and events simulated from LIRG (solid line). Taking into 
consideration the sum of all cosines above 0.99, the resulting value lies within 1.67 
a from the mean value obtained from the LIRG scenario, confirming quite a strong 
correlation between the real and simulated LIRG events. 

In figure 7 there are presented distributions of 47 events from AGASA and 
equivalent samples of events simulated from the map in figure 2 (top). We have 
estimated that the probability of getting a doublet or a triplet in the vicinity of Arp 
299 is about 0.5. We can see that samples of the 47 events simulated from LIRGs 
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Figure 7. Sky distribution comparison of the 47 AGASA showers with energy 
40-80 EeV (top left) and examples of three simulations from LIRG anisotropy 
map. Equatorial coordinates. 



anisotropy scenario, show good resemblance to the experimental data distribution. 
4. Astrophysical implications 

The strongest evidence for UHECRs origin may come from the correlation between 
directions of the AGASA triplet of events with energies 53.5, 55.0 and 77.6 EeV and 
an energetic astrophysical source in the local Universe. So far Arp 299 is the best 
candidate as a member of LIRGs, the brightest infrared source within 70 Mpc and 
a system of colliding galaxies showing intense, violent starburst activity at only 42 
Mpc away. It should be noted that Arp 299 appeared earlier as VV 118 in the list 
of candidates for the AGASA triplet presented by Takeda et al |^ but was not 
recognized as a colliding, energetic system, and, as a result, was not given enough 
attention. These authors point to another object VV 141 being a colliding system 
at z=0.02. However, as Arp 299 is at a distance two times smaller and fulfils the 
necessary criteria for CR acceleration, we think that it is this object which could be 
the most probable CR source. 

Recent studies on the nature of LIRGs suggest that compact strong radio emission 
may result from frequent multiple luminous radio supernovae pl[ . These objects are 
a poorly known class of supernovae with high nonthermal radio power indicating large 
kinetic energy input to accelerate particles p9| . In systems with intense starburst 
activity in extremely dense molecular gas and strong magnetic field environment, 
"nearly every supernova explosion results in a luminous radio SN with very high radio 



power" |21 . It is also worth to note that there are some observations suggesting a 



relationship between gamma ray bursts (GRB) and supernovae explosions |16| . Very 
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recently Weiler et al ||2^ have stated that GRB980425 and radio loud supernova 
SN1998bw are possibly related. Thus, it is not unreasonable to invoke here the 
intriguing hypothesis of the common origin of the two most energetic, mysterious 
phenomena in the Universe, UHECRs and GRBs. 

5. Conclusions 

We have shown that the available data from AGASA in the energy range 40 to 80 
EeV are not in contradiction with the expected anisotropy of CR produced in LIRGs. 
After applying the GMF, the LIRG hypothesis passes as well, but an apparent worse 
agreement with the data, suggests that the GMF model used by us may not be 
appropriate. So, we may hope that, if point sources existed they would be useful 
for determining the global structure of the GMF. At energies above 80 EeV both the 
isotropic and the LIRGs anisotropic distribution hypotheses are rejected. This might 
be explained by an existence of a UHECR population of a different origin above GZK 
cut-off. 

However, it seems that the SCvM test used here is not conclusive for distinguishing 
between the isotropy and the LIRGs distributions. Contrary to calculations done by 
other author, the eigenvector analysis has also appeared unable to do this for the 
existing sample sizes. Thus, we have considered a small scale clustering rather than 
a large scale anisotropy. Indeed, the analysis of multiplet probability gives good 
discrimination between isotropic and anisotropic scenarios (figure 5). The probability 
of the occurring of two doublets and one triplet (three doublets and one triplet) is 10 
(20) times higher for the LIRG hypothesis than that for isotropy . The analysis of the 
distributions of the angular distances between the data and simulated LIRG events 
also supports an existence of a correlation (figure 6). The strongest argument for our 
hypothesis is the observation of the triplet from the direction of Arp299, a system of 
colliding galaxies. We have obtained, via simulations, a high probability of doublets 
and triplets from this source, estimated to be about 0.5. 

Finally, with the prospect of new data to come from the Pierre Auger Observatory 
in a few years, the arrival direction distribution of UHECR should provide a much 
better clue to their origin. 
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Appendix A. Small angle scattering 

Let us consider the statistical process of the small angle scattering of a charged particle 
on N randomly oriented magnetic cells. Here we derive formulae for deviation angles 
and time delays on simple assumptions that statistical variables Sfik (denoting change 
of unit vector along direction of flight occurring in the fc-th cell) are: 1) independent 
and small, 2) cumulative change of direction is also small, 3) the number of cells N is 
large. As the particle energy changes along its trajectory, we have the relation 
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Figure Al. Schematic diagram of the particle trajectory through IGMF 




(^".) = (^^■) • I ^ I =<^^^V(^^j (A.1) 

Here En and denote particle energies in TV*'' and fc*'' cells and {•&%) is the mean 
square scattering angle in the last Nth cell, which is at the detector. We also denote 
cells by a continuous variable u = k/N instead of k (i.e. u = 1 indicates the detector). 
In this model the particle's trajectory consists of N segments of straight lines of equal 
length L = 1 /N. If the source of particles is placed in the centre of a sphere of radius 
R = I, then particle trajectories will end inside the sphere and the time delay is just 
the additional time the particle needs to reach the sphere. Looking back from the 
detector we miss the source by the deviation angle t] (figure Al). The formulae for the 
coordinates of the trajectory end: x,y, z (z-axis is the particle initial direction), time 
delay T and angle r\ are shown below. The cumulative change of direction is described 
by the angle 9^ as: nl - rTo = X^*^^ fei = Ym=\ '^i = ^k- 

1 ^ 1 ^ 

fc=l fc=l 



z = 



fc=i fe=i 



where Az = - ■ Y^h^lk + A) (A.3) 
fe=i 
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We can calculate statistical moments for the above random variables. Here are some 
of them: 

{rf') = N ■ i^r,') ■ /12; iff) {ff? = ivY 

(r> = y • {^n') ■ (/2 - /ii) (A.4) 
The constants /n, /12 and /2 are defined by the following integrals: 

u^A{u)du; 



hi = / (1 - u)^A{u)du; /12 = / ' 



/2 = / {l~u) A{u)du (A.5) 



where 







Aiu)^\^] (A.6) 



is the particle energy relation from source to detector assuming continuous energy 
losses. It should be noted that allowing for fluctuations in the energy losses would 
give larger deflection angles. The values of the derived moments have been positively 
checked in simulations with various functions A{u). 

Appendix B. Probability of multiplets. 

Here we derive analytical formulae for average numbers of doublet, triplet and 
quadruplet events on the assumption that the expected flux is known for any direction 
and single showers occur independently of each other. The last assumption means that 
a particular number of multiplets of a given kind (e.g. that of doublets) undergoes 
Poisson distribution with the expected value equal to the average of the expected 
(local) values over the entire map. Let us denote by X{S, a) {6 - declination, a - RA) 
the expected angular density of showers, i.e. the expected number of showers per unit 
solid angle within the measurement time (A — const for isotropy). The actual number 
of registered showers depends, of course, on the efficiency ri{d, a) of the particular air 
shower array. Assuming that only showers with zenith angles Q smaller than some 
Qz are considered, we define here that j]{d,a) — {f{t)cos{Q{t))) where f{t) = 1 while 
Q{5,a,t) < 8z and f{t) ~ otherwise. The brackets mean time average. It can be 
derived that: 

77(5,0;) = [sin(iy9) sin(5)Q;2 + cos((y9) cos((5) sin(Q!^)] /tt (B.l) 
where the angle is defined by: 

cos(e^) - sin(y) sin((5) 
cos a^) = — — B.2 

cos{ip) cos(d) 

This formula agrees well with the experimental efficiency of the AGASA experiment 
[ p5[ . The expected angular density p{6, a) of registered showers equals r]{S, a) ■ A((5, a) . 
If the mean number of showers detected is N, then we have the normalizing relation 
N = J p{5, a)dn and the average flux A is determined: 

N 

A = r , (B.3) 
J rj[d,a)d\l 
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We assume that the number of registered showers in an experiment is large enough to 
adopt it as the average number of showers which would have been registered by the 
detector in many such runs. The criterion for classifying two showers as a doublet is 
that the angle between them should be smaller than some given small value R. As 
the expected number of showers within the small angle R around a given direction 
is ^ = p{ttR^) (assuming constant p within the circle), the local detection density of 
doublets can be written as p2 — ^ppe^f^ . The factor i stands for the fact that both 
members of a doublet count as one doublet. Of course, the local detection density of 
singles, where there is no companion within the angle R from the shower considered 
is pi = pe~^^. 

The situation is more complex for the case of a triplet. Two criteria may be 
distinguished: 

1st - another two showers should be closer than R to the shower considered (in the 

centre) - criterion applied by AGASA. 
2nd - three showers should be within a circle with diameter R. 

The 1st criterion is illustrated in the figure Bl (left). In the centre of the right circle 
there is the shower denoted as A and within the radius R there are another two 
showers, out of which the shower denoted as B is in a distance x from the first one, 
in the centre of the left circle. The third shower may be either on the area where the 
two circles overlap (Ca) or outside it (Cb). In the first case around all three showers 
a triplet will be found (the same) . In the second case the triplet is found only around 
shower A. 

Let S(x) denote the ratio of overlapping area on the figure to the area of the circle. 
We have the following relations: 

, , i?2r2a- sin(2a)l , ^ 

Si^) = — p2 B.4 

where cos(q:) = x/{2R). The average value of the ratio S{x) equals: 

{S{x)) = S{x)^dx = 1 - ^ (B.5) 

Thus, the local detection density of triplets under 1st criterion can be written as : 



P3 = P 



\{S{x)) + (l-{S{x))) 



^e- 
2! 
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The mean numbers of doublets N2 and triplets detected from all directions are 
found by integrating their detection densities over the whole solid angle. 
We get: 

iV2 = P2{S,a)dn = ^(Ti-i?^) j {Xv^e'^dn w ^{-kR^)X^ j rfdO. 
Ns^ J p3{S,a)dQ = i(7ri?2)2 (^i + J {Xr,fe->^dn 

^^-inRrx'(^l + ^^Jv'dn (B.7) 

where we have assumed that A = const. Similar consideration as those for detection of 
triplets under the 1st criterion lead us to the corresponding formulae for quadruplets 
under the same criterion: 



Pi = P 



1 „ 1 „ „ „ N 1 /i^ 



3!" 



j P4{S, a)dn = ^{nR^ fX'^ ■ 0.66 • j jy^e'^dl^ (B.. 



where S34 is the probability that three showers randomly distributed within a circle 
with a forth shower in the centre will constitute the same quadruplet if any shower 
out of the three is chosen as the central one. The probability that only one shower out 
of the three, chosen as the central one, makes the same quadruplet possible is denoted 
514. It has been found by integration that 5*34 = 0.274 and 514 = 0.270. 

The formulae for detection densities P2,P3, P4 derived above allow for such shower 
configuration, where the shower B in the figure Bl (left) being a member of a doublet 
or a triplet is at the same time a member of a triplet or a quadruplet, respectively. 
After excluding such situations we get for doublets the following formula: 

p^ = p,-p.3{l-{S{x)})(^e-'^ (B.9) 

Figure Bl (right) illustrates the situation for classifying three showers as a triplet for 
the 2nd criterion, i.e. three showers have to be within a circle with diameter R. The 
two circles there have diameter 2r = R. Two showers from directions A and B are at 
the angular distance AB < R. To find the area where the third shower should be we 
draw the circle of radius R/2 through the point A and rotate it around point A until 
point B is inside the circle, so that points E goes to point E'. Similarly, we rotate the 
circle around point B until point A is inside the circle (point D goes to D'). To form 
a triplet the third shower must be located within area E'EDD'. 

Let S2{x) be the ratio of the framed area to nR'^ for a given value AB= x. It 
follows that: 

1 X 

5*2(2;) = |"6q: + TT — sin(2Q!)l — where cos(q:) = — (B.IO) 
'■ 47r R 

Taking into account the distribution of distances AB we can calculate the average 

value of the 52 (x): 

{S2{x)}=2£s2{x)^dx = ^ (B.ll) 
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Further considerations for the local detection density of triplets under the 2nd criterion 
are only for the case /i <C 1. So, in a similar way, we get that: 

p, = pl{S,ix))^^^p^,' (B.12) 

Numbers of observed doublets, triplets, quadruplets etc undergo Poisson distributions 
with the corresponding mean (expected) values 7V2,A^3, {N2,N^), etc as calculated 
above. The formulae for the mean numbers of the multiplets have been positively 
checked in simulations. For the AGASA experiment the following values apply: 

= 45°;^ = 35°47';i? = 2.5° 
We find that : 

/ r]dn = 7r(l - cos^e^) = 1.57 sr; / r]^dn = 0.38 sr; / 7]^dn = 0.097 sr; / r/'^dQ = 
0.0255 sr 

For the expected value of showers N we use the actual number (registered in the energy 
range 40 to 80 EeV) - 47 showers and obtained that: 

N2 = 1.018; N3 = 0.0283; = 0.954; N'^ = 0.0268 under 1st criterion; N3 = 0.0087 
under 2nd criterion; A'4 ~ 4.57 x 10^"* for quadruplets under 1st criterion. 

Thus the probability of obtaining two doublets or more and one triplet or more, 
for the isotropic sky without point sources equals [l — ea;p(— iVj) — N2 ■ ea;p(— iVj)] • 
[1 - expi-Ns)] = 0.247 • 0.0279 = 6.9 x 10"^. We see that it is mainly the triplet 
event that makes the isotropy hypothesis very unlikely. 
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